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Abstract 

We investigate the mechanism that leads to systematic deviations in cluster Monte 
Carlo simulations when correlated pseudo-random numbers are used. We present a 
simple model, which enables an analysis of the effects due to correlations in sev- 
eral types of pseudo-random-number sequences. This model provides qualitative 
understanding of the bias mechanism in a class of cluster Monte Carlo algorithms. 
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1 Introduction 



The Monte Carlo method to obtain statistical averages for a model is widely 
used when exact calculations are not available. The error, expected from a 
standard statistical analysis, is proportional to the inverse square root of the 
number of randomly chosen samples. One may find biased results if the samples 
are not representative. Upper bounds to the error are not known in general. 
Standard statistical analysis of simulation data suffices to obtain the errors 
of a random statistical nature. Statistically independent random numbers do 
not introduce a bias in Monte Carlo simulations. However, it is not possible 
to achieve independence with a deterministic recipe, as is commonly used in 
Monte Carlo calculations[l]. Although one generates seemingly irregular num- 
bers the underlying production rule induces relations or correlations between 
the pseudo-random numbers. The production rule itself already constitutes a 
correlation between the generated pseudo-random numbers. The correlations 
may cause a bias in the simulation results that is difficult to assess. 

Substitution of the production rule in itself leads to further correlations involv- 
ing more random numbers. By repeated substitutions one finds a hierarchy of 
correlations. Such hierarchies were studied by Compagner[2,3]. Due to the mul- 
titude of possible correlations it is not feasible to study all of them. One has to 
analyse which correlations are detrimental to the intended use of the pseudo- 
random numbers. Compagner notes, that correlations between few numbers or 
closely spaced numbers arc most important. The effect of such correlations has 
actually been observed in Metropolis type Monte Carlo simulations of Ising 
models [4,5]. 

In this work, we investigate the consequences of correlations between pseudo- 
random numbers in an example of a simple random-walk model. An analysis of 
the effect of the correlations present when a given generator is used may guide 
us in the choice of an appropriate production rule. In Section 2 we introduce 
our stochastic walk model. In Section 3 we discuss the class of generators we 
used. Deviations caused by the correlations are described in Section 4 for shift 
register generators and in Section 5 for the lagged Fibonacci random number 
generator. In the last Section we give a discussion and conclusion. 



2 The walk model 

The Ising simulations using cluster updates with the Wolff [6] method seem 
to be very sensitive to deficiencies in a number of pseudo-random-number 
generators [7,8]. The crucial element of such simulations is the formation of 
clusters. These are formed by joining spins which are connected by 'active' 
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bonds. In ferromagnetic models bonds can only be active, if the spins on 
either side of the bond are in the same state. The probability, that a bond is 
active, still depends on the coupling between the spins. 



In order to study the random-number generator induced bias we replace the 
Wolff cluster formation process by a simpler one. First we place the Ising spins 
on a one-dimensional lattice. Second the cluster is grown only in one direction 
(or, equivalently, the cluster formation always begins at an open boundary of 
the spin chain). Third the cluster is grown in a configuration of parallel spins 
only. This third simplification is less far reaching than it may seem. It does 
not modify the Ising statistics of the system. The probability that Ising spins, 

coupled with a strength K, are parallel is (l + e~'^^^ . If they are parallel the 

random-cluster model yields a probability (1 — e~^^) that they are connected 

by an active bond, i.e. that the two spins belong to the same cluster. Thus, in 
Woffi simulations of the one-dimensional Ising model the probability that each 
next spin is included in the cluster, equals //. = tanh(i^): a constant depending 
on K only, just as in our simplified model. 

We may interpret the cluster formation in the simplified model as a directed 
one-dimensional random walk. A walk starts on one site (see Fig. 1). At dis- 
crete times n the walker steps to the neighbouring site in a fixed direction (to 
the right in Fig. 1) with probability /j,. Thus the walker cannot recur to a site 
once visited. If the step is not made, a new walk is initialized. The decision, 
whether to start a new walk or not, does (ideally) not depend on previous de- 
cisions. The probability to start a new walk at a certain time equals u — 1 — //, 
independent of the length of the walk. The probability P{n) that this process 
generates a walk covering precisely n sites, satisfies 



One can formulate this random walk process in terms of the following two 
algorithms. We distinguish two possible ways to generate the decision whether 
to step to the next site or to initialize a new walk . The algorithms differ by 
the underlined statements. 



P{n) ^ - fx) 




and the expectation value of the number of sites visited is 




(2) 
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Algorithm P 

Initialize walk statistics {ul} 

generate a random number r 
(with < r < 1) 

for walk := 1 to total number do 
{ walk length L:— 1 
while > r do 
{ L:=L + 1 

generate a random number r } 
enddo 

riL := riL + I } 
enddo 

output {ul} 



Algorithm N 

Initialize walk statistics {til} 

generate a random number r 
(with < r < 1) 

for walk := 1 to total number do 
{ walk length L:— 1 
while 1 // < r do 
{ L:=L + 1 

generate a random number r } 
enddo 

ul -.^ Ul + I } 
enddo 

output {ul} 



Using perfect random numbers, independent and uniformly distributed be- 
tween and 1, both algorithms would yield consistent and unbiased averages. 
However, in real simulations, the random numbers arc correlated, as a con- 
sequence of the deterministic production rule. This causes the distribution 
of n-tuples of consecutive pseudo-random numbers to be non-uniform in the 
unit hyper-cube especially for larger dimensionality n. A certain degree of non- 
uniformity is acceptable for the numbers used in the directed random walk al- 
gorithm. It is sufficient, that n-tuples with components larger or smaller than 
the quantity used for the decision (/x in Algorithm P and 1 — in Algorithm 
N) are present in the right proportions. But even this weaker requirement is 
not satisfied in practice. This affects the expectation values of the simulation 
results. These deviations are discussed in Sections 4 and 5. 



3 Random generators 



In a numerical Monte Carlo calculation as described in Section 2 one has to 
decide between stepping with probability ^ and initialization of a new walk 
with probability u = 1 — fi. One needs a random sequence of two decision 
outputs, where one occurs with the step probability and the other with 
the initialization probability 1 — To this purpose we use pseudo-random 
numbers r with < r < 1 and compare them to the relevant probability. 

We generate the pseudo-random numbers with a number of rules. One is the 
Generahzed Feedback Shift Register (GFSR) method [9,10]. A sequence of 
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pseudo-random integers r{i), represented by their binary expansion (usually 
32 bits), is generated by the rule r{i) = r{i—p)®r{i — qi)®- ■ ■(Br{i — qj), where 
the exclusive-or operation © is applied bit-wise and the feedback positions 
(lags) are ordered according to p > qi > q2 > • ■ • > qj. li j = I, we denote 
the feedback lag by q. The leading bits (as well as all other bits) separately 
form a sequence, that is generated by a feedback shift register. The maximum 
length equals 2^ — 1. Production rules that generate maximum-length bit- 
sequences, are characterized by primitive polynomials [11]. Tables of primitive 
polynomials are listed in [12]. The most widely used random number generator 
based on the GFSR method is R250 [10] with p = 250, q = 103 and r{i) = 
r(i- 250)©r(i- 103). 

Where primitive trinomials are known, one can obtain other primitive polyno- 
mials by decimation procedures [11]. For example, using every third number 
from a sequence generated with a rule derived from a trinomial is equiva- 
lent to using a rule specified by a pentanomial, which adds two terms to the 
primitive trinomial. For a rule based on primitive polynomials of degree p 
only sequences up to p successive bits or numbers are independent. In order 
to prevent problems one should avoid unwanted initializations. For instance, 
if the k-th bit of each initialized integer is zero, the same will hold for the 
subsequently generated pseudo-random integers. 

Furthermore we shall consider the lagged Fibonacci method: r{i) — r{i — p) + 
r{i — q), where the addition is understood modulo 2', where I is the number of 
bits in a computer word. The maximum sequence length[13,14] is 2'~^(2*' — 1). 



4 Examples of the correlation effect: Shift Registers. 

An important factor contributing to the bias for a given random-number gen- 
erator is associated with the fact that, whenever a new walk is initialized, 
the preceding decision was not to visit another site. The last random number 
rp(0) used in the previous walk in Algorithm P thus obeys rp(0) > (or obeys 
?^Af(0) < 1 — jj, ~ u in Algorithm N). If the new walk visits another site, the 
next random number rp(l) < n (orr7v(l) > ly)- If the computer words consist 
of / bits, only probabilities that are a multiple of 2~^ are faithfully represented. 
We suppose that the first p numbers of a generalized shift-register sequence 
may be considered as independent (as if produced by an ideal random number 
generator), i.e. the decisions with respect to the first p — 1 steps occur with 
probability // independent of the history. The probability of a walk with less 
than p steps is that given by Eq. (1). 

However, deviations will occur for step number p. The p-th number in the 
sequence depends on the bits in r(0) and those in the integer r{p — q) on the 
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feedback position. This causes the initiahzation probabihty at this point to 
differ from the value v it should have. We denote the actual probability by 
V* . This probability will depend on the chosen random generator (GFSR or 
lagged Fibonacci, etc.). In this Section we will focus attention on shift-register 
random-number generators. The case of the lagged Fibonacci recipe is treated 
in the next Section. 

First, we analyze some simple examples. \i \i = \ just one random bit h 
is needed to generate the decision, where both outcomes occur with equal 
probability. The random bit h denotes the leading bit of the pseudo-random 
number r. The completion of the preceding walk implies, except for the first 
walk, that 6(0) = 1 when Algorithm P is used. In the case that a decision is to 
be made on step p, the completion of the preceding steps means that h{i) = 
for < i < p. Thus, irrespective of the number j -|- 1 of feedback positions 
P,Qi,Q2, ■■■,Qj, b{p) = 1 so that the walk ends, and a new one is initialized. 
The actual initialization probability is denoted i'*{p) and satisfies i^p{p) — 1 
instead of the desired u = 1 — fi. The maximum walk length is p and occurs 
with a probability of twice that expected from Eq. (1). 

In the case that a decision has to be made on step p of Algorithm N, 6(0) = 
and b{i) = 1 for < i < p. If the number of additional feedback bits is odd, 
which is the case for maximum-length sequences, b{p) = 1 and the walker 
always proceeds to site p. The actual initialization probability thus satisfies 
^jv(p) = 0- However, one has b{p + 1) = 0, so that i^^iP + 1) = 1- Thus all 
walks with a length exceeding p — 1 have a length p + 1. Length p occurs with 

probability 0, length p + 1 with probability P*{p+1) — (^^Y . 

These results can be generalized for all n = 2"™" and positive integer m. The 
walk lengths will be smaller than p + 1 in Algorithm P, as the leading m bits 
of r{p) will be the same as those of r(0). For Algorithm N the leading m bits 
of r{p + 1) have to be zero. Lengths larger than p -|- 1 do not occur. 

If /i > |, much longer walks may occur, if the GFSR-generator is used. As 
an example of this regime we take /i. = |. One needs three bits to make a 
decision. The possible values are represented by the integers f = [8rJ where 
the brackets denote the integer part. 

In the case of Algorithm P a decision on step p requires f(0) > 5 and r{i) < 5 
for < i < p. We assume that all admissible numbers for f{i) with i < p 
occur with equal probability, which strictly holds only for the first walk. For 
sequences generated with a production rule derived from a primitive trinomial, 
f(p) equals 4, 5, 6 and 7 with probability |, 1, 2 and 3 with probability 
and with probability 0. This non-uniform distribution of f leads to a 
probability ^'p(p) = | for the walk to end at step p. This is different from the 
desired probabihty u — ^. The probabihty of walks of length p thus equals 
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,p-l. 



For r(i) with 



p <i < < 



p + q if g < f 



2p-qii q> ^ 



(3) 



the numbers 4, 5, 6 and 7 occur with probabihty ^, the numbers 1, 2 and 3 
with probabihty ^ and the number with probabihty |. Then ^'p(i) = ^. 

Chains of length p+ 1 occur with probability P*{p+ 1) = (|)^ | ■ ^ provided 
q 7^ 1, p — 1. Similar non- uniformities of the distribution lead to i^p{p + q) = | 
or iy*p{2p-q) = |. 

For a production rule derived from a primitive pentanomial with lags p > qi > 
<l2 > Qs the actual initialization probability i^p{p) = and i'p{i) = ||| for i 
larger than p and less than the minimum of p + 53, p + ^2 — qs, P + qi — q2 and 
2p — ?i- The deviations are less than for trinomials. 

Whenever a decision on step p is made in the case of Algorithm A^, the pseudo- 
random numbers satisfy f(0) < 2 and r{i) > 2 for < i < p. Assuming that 
f(i) with i < p occurs with equal probability for all possibilities consistent 
with a walk length larger than p — 1, f(p) equals 4, 5, 6 and 7 with probability 
|, 1, 2 and 3 with probability ^ and with probability in the case of a 
primitive trinomial. So ^'^(p) = ^5. The probability of walks of length p thus 

equals P*{p) — (|)^ j^- For f(i) with i obeying inequality (3) the values 4, 
5, 6 and 7 occur with probability the values 1. 2 and 3 with probability ^ 
and the value with probability i. Thus ^'iV(i) = Chains of length p + 1 

occur with probability P*(p + 1) = (^^Y ^ • ^, if g 7^ l,p — 1. In a similar 
way i/^(p + g) = § or z/^(2p = 

In the examples a resonance in versus i is found, when the value of r(0) 
affects the value of r{i). The first one encountered is a reflection of the produc- 
tion rule itself; for trinomials a relation between lags 0, q and p. The next one 
is a four-point correlation resulting from the interference of the production 
rule and a shifted version of the production rule. The correlation between lags 
0, q and p and that between lags q,2q and p-|- g leads to a correlation between 
lags 0, 2q, p and p + q. Similarly the correlation between lags 0, q and p and 
that between p — q, p and 2p — q leads to a correlation between 0, p — q, q 
and 2p — q. The deviation of the initialization probability has a different sign 
for Algorithm P and Algorithm A^, if the new correlation was between an odd 
number of lags. It has the same sign for the four-point correlation. 
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If a walk is longer than p in the case of Algorithm P for a general real step 
probability ^ > |, then rp(0) > ^ and rp{i) < fi for < i < p. As a 
consequence of these inequalities, the production rule generates r(p) with a 
non-uniform probability distribution. For a GFSR-rule derived from a trino- 
mial r(p) — r(0) ® r{p — q). The exclusive-or operation with a fixed operand 
is a permutation of the numbers smaller than 1. If both r(0) and r{p — q) are 
larger than // the result is less than |. This implies that all possibilities for 
r(p) > 11 are realized for the reduced set of numbers r{p — q) < ji. So the prob- 
ability that the walk ends, when it has visited p sites, is not equal tov — 1 — /x, 
but Vpip) — ^. By similar arguments for i obeying (3) i'p{i) — {'ijJi — ji"^. 



For Algorithm the number rjv(O) < v and rN{i) > v. \i v = 2^™ with m 
a positive integer, then rAr(O) has m leading bits equal to zero, so rAr(p) has 
the leading m bits the same as r^ip — q) and ^'Xr(p) = 0. A similar result 
is found for all lengths for which a three-point correlation involving rAr(O) is 
induced by the production rule. In particular in the case of a rule derived 
from a primitive trinomial this holds for all lengths 2^p with k a non-negative 
integer, because a decimation by 2 leads to a sequence generated by the same 
rule [11]. For general v a more complicated argument leads to 

and for i obeying (3) 

with 



no = riog2(i^2' + 1)1 - 1. 



where \x\ is the smallest integer > x. 

The deviation of the actual initialization probability u* from the ideal value 
V causes the probability of a walk to visit n sites P*{n) to deviate from the 
value for uncorrelated numbers (1). We define 

SP{n) ^ ^ - 1, P*{n) = n(l - ^*«) (4) 



For i < pwe expect no deviation, so SP — 0. Substitution in (4) of the values 
found for u* leads to 5Pp(p) — ^ for Algorithm P. With i obeying inequahty 
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(3) 



(2/^-1)2 /3/i2-3/i+l 



) 



i—p—l 



5Pp{i) 



- 1 



To give an impression of the deviations for lengths beyond the first new low- 
order correlation we use numerical calculations. The results of computer simu- 
lations of 10^ walks are presented in Fig. 2 for the GFSR with {p, q) = (89, 38), 
Algorithm P and /i = 31/32. We mention two reasons to choose these partic- 
ular values of p and q. First, simulations with the desired accuracy typically 
need more than 10^^ random numbers. Therefore, the length p of the shift reg- 
ister should be greater than log2 10^^ 40. Second, we do not want to choose 
p much larger than necessary because the observability of the bias decreases 
with p. 

The error bars were computed using 100 bins of 10^ walks each. The predicted 
resonances are easily seen for the walk lengths p, p + q, 2p — q, 3p — 2q, etc., 
which are linear combinations of the feedback positions in the production 
rule. The result for walk length 2p is not shown in the Figure because of its 
large deviation of 5P{2p) = 0.109(3). The value of the calculated deviation 
SP{p) = 0.0328(7) is in good agreement with the value ^ ?s 0.0323. 

The results of the calculations with the same parameters using Algorithm N 
are shown in Fig. 3. The deviations at lengths equal to 2^p (A;=0,l,2...) are 
equal to —1 and are not shown in the Figure. The deviation 5P{p + 1) ~ 
0.0651(8) is in agreement with the value (1 — /x^)///^ ^ 0.0656. A comparison 
of Fig. 2 and Fig. 3 shows, that the main resonances have different signs for 
the two algorithms. This is in qualitative agreement with the results of cluster 
simulations of the two-dimensional Ising model by Shchur and Blotc [8] who 
indeed found deviations in thermodynamical quantities with opposite signs 
for comparison strategies P and N. 

In order to give an example of the deviations of random walk statistics for 
a more complicated production rule we use the decimated sequence with 
(P)9i)?2,?3) = (89,72,55,38), Algorithm P and = || in a random walk sim- 
ulation. Because the feedback positions for this sequence are equally spaced, 
the four numbers with lags 106, 38, 17 and are correlated. This four-point 
correlation may be more important in applications than the five-point correla- 
tion implied by the production rule. This is the case for the numerical results of 
the random walk statistics (Fig. 4). The deviation at the shift-register length 
SP{p) — 1/31^ is small compared to the error bars. The deviations for the 
two four-point correlations at walk lengths 106 and 212 are of a similar mag- 
nitude as the deviations for four-point correlations for the trinomial (Fig. 2). 
The number of four-point correlations is smaller than for the trinomial and 
no three-point correlations are found. Deviations tend to be weaker for higher 
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order correlations. This agrees with resuhs of Wolff simulations of Ising mod- 
els[8], where the deviations were smaller for decimated sequences. 



5 Correlation effects due to the lagged Fibonacci recipe. 

For lagged Fibonacci generators the relation between the leading bits is less 
simple. The distribution of numbers is not completely uniform for lagged Fi- 
bonacci sequences. Deviations of uniformity have been analysed in [15]. We 
assume that all sequences r(0) . . . r(p — 1) are equally probable. If the walk has 
a certain length the previous numbers r{i) obey the appropriate inequality. 

Wc first consider /U = | again. For Algorithm P, r(0) > | and r{i) < ^ for 
< i < p, in the case of a decision on step p. These inequalities cause the 
production rule to generate r(p) with a non- uniform probability distribution. 
The distribution can be calculated as the convolution of the distributions 
for the feedback lags. For numbers of infinite precision it has the symmetry 
property P(r) = P(l — r) in the case of a rule derived from a primitive 
trinomial. The finite word length / leads to corrections of order 2~K We neglect 
these corrections because they are small compared to the statistical errors in 
our simulations. 

The actual initialization probability conserves its ideal value i^p{p) = | by 
the symmetry property. The same symmetry holds for the distribution of r(i) 
with i obeying inequality (3). Thus ^'p(i) = ^. For the next random number 
the probability distribution does no longer obey P{r) = P{1 — r). The actual 
initialization probability z/p(p + g) = | or i'p{2p — q) — ^- The same results 
are obtained for algorithm N. 

As in the case of the GFSR rule it is possible to calculate the values of the 
actual initialization probability for arbitrary /x. For a lagged Fibonacci rule 
derived from a trinomial the actual initialization probability becomes i^pip) — 
For i obeying inequality (3) the actual probability ^'p(i) = (3/x — l)i//2/x^. 
For strategy N we get the same actual probabilities. The expressions for u* 
grow progressively more complicated for larger walk lengths. We therefore 
refrain from showing them. 

Numerically calculated deviations of the resulting probability P*{n) of a walk 
of length n are shown in Fig. 5 for the lagged Fibonacci rule with feedback po- 
sitions (p, q) — (89, 38). Error bars were computed using the same parameters 
as in the previous Section. Because of the number of bits / = 30 used in the 
random numbers the finite-word-length corrections are quite small. As those 
differences are small compared to the statistical errors, the results for Algo- 
rithms P and N are equal within error bars for all walk lengths, unlike in the 
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case of shift registers. The value of the deviation at walk length p (not shown 
in the Figure) is equal to SP{p) = —0.4841(4) in good agreement with the 
value of ^ — 1 a; —0.4839. The deviations have the same order of magnitude 
for both recipes of random number generation we considered in this paper. 
This does not support the claim that the lagged Fibonacci method performs 
better than the GFSR-method[14]. 



6 Conclusions 

Correlation between random numbers can influence the results of a Monte 
Carlo calculation of a simple random walk model. No deviations occur in 
the distribution of walk lengths shorter than the magnitude of the largest 
feedback position. The walk length statistics is affected for larger lengths. In 
some cases the difference in results using comparison strategy A'^ and P gives 
an indication of the magnitude of the bias. The magnitude of the deviations 
tends to be larger for generators derived from primitive trinomials than for 
primitive pentanomials. Similar effects occur in the cluster-size distribution in 
a cluster simulation of the Ising model[16]. 

The magnitude of the deviations i/ — i/* depends only on the comparison strat- 
egy, and on the value fi. It does not depend on the particular values of feed- 
back positions. In the case of other feedback positions, for example R250 the 
Kirkpatrick-StoU Random Number Generator, one should accordingly rescale 
the horizontal axes of Figures 2-4 by 250/89. However, continuation of the 
walks to larger n leads to increased scatter so that the effects become less 
prominent. 

These results are relevant for the Wolff cluster simulation of spin models in 
more than one dimension. In such simulations, each spin in the cluster may 
have more than one neighbour that has to be included in the cluster. If so, the 
addresses of these neighbouring spins are temporarily stored, e.g. in a memory 
called 'stack' ( for a hardware implementation see [17]). For each entry in the 
stack, it has to be checked whether further neighbours have to be added to 
the cluster, and thus whether further additions to the stack memory have to 
be made. Spin addresses that have thus been processed, are removed from 
the stack. The number of addresses in the stack is a fluctuating variable, and 
the cluster is completed when the stack is empty. Typically the stack memory 
contains more than one spin address, in which case one random number is not 
sufficient to end the cluster formation process. However, it is obvious that the 
completion of a cluster is strongly correlated with the values of the preceding 
pseudo-random numbers. 

The bias in the random numbers preceding the construction of a new cluster 
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will lead to a significant correlation between n — 1 pseudo-random numbers 
in the case of an n-point production rule, although the correlation is weaker 
than the correlation between n numbers imposed by the production rule. For 
a 3-point rule this 2-point correlation combined with the bias in the random 
numbers allowing a further growth of the cluster causes a bias for the p-th 
pseudo-random number used in the construction of a cluster. It is known[8] 
that the bias in the simulation results (for three-point production rules) be- 
comes largest, when p is equal to the average cluster size. 

For higher n the bias-producing mechanism is less simple and it seems plausible 
that the bias in the p-th pseudo-random number will be weaker. Thus one 
expects that the deviation in the statistics of cluster sizes will be stronger, 
when less pseudo-random numbers are involved in the production rule. Indeed 
for the case of shift-registers with a n-point production rule, the bias in the 
thermodynamics in a Wolff simulation increases strongly with decreasing n [8], 
where rules with higher values of n can be generated either by decimation or by 
combining the numbers generated by two or more rules with an exclusive-or. 
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Fig. 1. Graphical representation of the walk generation. A step is made with prob- 
ability fj,, a new walk initiated with probability 1 — /i. The random number r{k) is 
used to decide whether the step is made from k to k + 1. 
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Fig. 2. Deviation SP of the probability of a walk of length n from the value for 
uncorrelatcd random numbers versus walk length for a Monte Carlo calculation 
using Algorithm P, ^ = |^ and GFSR-rule (p,q) = (89,38). Some of the resonances 
are labeled by linear combinations of feedback positions. The deviation at length 
2p is equal to SP{2p) = 0.109(3) and outside the scale of the figure. 
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Fig. 3. Deviation 5P of the probability of a walk of length n from the value for 
uncorrelated random numbers versus walk length for a Monte Carlo calculation 
using Algorithm N , n = ^ and GFSR-rule (89,38). The deviations are -1 for length 
p and 2p. 
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Fig. 4. Deviation 5P of the probability of a walk of length n from the value for 
uncorrelated random numbers versus walk length for a Monte Carlo calculation 
using Algorithm A^, /x = § and GFSR-rule (89,72,55,38). 
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Fig. 5. Deviation 5P of the probability of a walk of length n from the value for 
uncorrelated random numbers versus walk length for a Monte Carlo calculation 
using Algorithm P, /x = |i and lagged Fibonacci rule (89,38). The deviation at 
length p = 89 equals -0.4841(4) and is not shown in the figure. 
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